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ABSTRACT 

We construct a two-parameter family of models for self-collimated, magnetized outflows 
from accretion disks. As in previous magnetocentrifugal wind solutions, a flow at zero 
initial poloidal speed leaves the surface of a disk in Kepler rotation about a central star, 
and is accelerated and redirected toward the pole by rotating, helical magnetic fields which 
thread the disk. At large distances from the disk, the flow streamlines asymptote to wrap 
around the surfaces of nested cylinders, with velocity v and magnetic field B directed in 
the axial (z) and toroidal (ip) directions. In the asymptotic regime, the velocity secularly 
decreases with cylindrical radius R from the inside to the outside of the flow because 
successive streamlines originate in the circumstellar disk in successively shallower portions 
of the stellar potential. In contrast to previous disk wind modeling, we have explicitly 
implemented the cylindrical asymptotic boundary condition to examine the consequences 
for flow dynamics. The present solutions are developed within the context of r- self-similar 
flows, such that v, the density p, and B scale with spherical radius r as v oc r~^/^, p oc r~'^, 
and B oc r^^^^''^^'^; q must be smaller that unity in order to achieve cylindrical collimation. 
We self-consistently obtain the shapes of magnetic field lines and the ^-dependence of all 
flow quantities. The solutions are characterized by q together with the ratios Ra/Ri and 
Rq/Ri, where for a given streamline Rq is the radius of its footpoint in the disk, Ra is the 
cylindrical radius where the flow makes an Alfven transition, and Ri is its final asymptotic 
cylindrical radius. For given q and Rq/Ri, Ra/Ri must be found as an eigenvalue such that 
the Alfven transition is made smoothly. In the solutions we have found, the asymptotic 
poloidal speed Vz on any streamline is typically just a few tenths of the Kepler speed VLRq 
at the corresponding disk footpoint, while the asymptotic rotation speed may be a few 
tenths to several tenths of VLRq. The asymptotic toroidal Alfven speed VA,ip = B^/ y/inp is, 
however, a few times QRq; thus the outflows remain magnetically dominated, never making 
a fast-MHD transition. We discuss the implications of these models for interpretations of 
observed optical jets and molecular outflows from young stellar systems. We suggest that 
the difficulty of achieving strong collimation in vector velocity simultaneously with a final 
speed comparable to QRq argues against isolated jets and in favor of models with broader 
winds. 
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1. Introduction 



1.1. Observational Context 



Energetic, coUimated jets and outflows are produced by accreting systems in a wide 
range of astronomical environments, from young stars forming in cold molecular gas clouds, 
to white dwarfs, neutron stars, and black holes in evolved binaries, to supermassive black 
holes within active galactic nuclei (e.g. p^ivio (1997)| ). Radiation pressure, thermal pressure, 
and magnetic stresses may all play a role in accelerating and coUimating such flows. For 
the case of low-mass pre-main-sequence (PMS) stars, the high mass and momentum losses 
observed ( |Lada (1985)]JEdwards, Kay, fc IVlundt (1993)| , [b'ukui et al (T993]i |Bachiller (1996)]) 
have led theorists to believe that magnetic forces are the essential ingredient in driving 
these winds (e.g. [Konigl fc Ruden (1993)| , |Shu (199^ ). The densities and velocities involved 
place such flows from young stars in the regime of nonrelativistic magnetohydrodynamics 
(MHD). The prevalence of jets and outflows from pre-main-sequence stars argues that 
they are an inevitable byproduct of star formation; outflows may in fact help determine 
the eventual mass of the star that forms ( phu, Adams, fc Lizano (1987)|) , and also have a 
profound influence on the dynamical evolution of the parent cloud ( |Norman fc Silk (1980)| ). 
Thus, understanding the physics of cold, nonrelativistic MHD winds from accreting systems 
is crucial for modeling of star formation, and also informative for studies of other systems 
where magnetic flelds are probably important in driving and conflning winds - but where 
electromagnetic, thermal, radiative, and relativistic effects may complicate the dynamics 
(e.g. IBegelman, Blandford, fc Rees (1984)^ |Blandford {im^ . 



A topic of great theoretical and observational interest in studies of PMS stars is the 
dynamical connection between the fast, ionized jets seen near the presumed polar axis of the 
wind, and the slower and less-collimated molecular outflows that surround them. Molecular 
outflows are sometimes "jetlike" in appearance ( [Bachiller fc Gomez- Gonzalez (1992)[ ) but 
often more poorly conflned spatially, while still having strongly directed bipolar momenta 
(e.g. [Lada fc h\ch. (199B}| ). Large outflow masses and red/blue lobe asymmetries argue that 
the bulk of the material in molecular outflows is swept up from the surrounding cloud rather 
than itself comprising a wind, but there is still no deflnitive model for the process that 
imparts momentum to the cloud. The observed line-of-sight-velocity/plane-of-sky-position 
distributions of outflow material must reflect the combined distributions of density in the 
ambient medium and forces or momentum fluxes that drive the molecular flow. Hence, 
although mass and momentum maps and line proflles from outflows are not directly 
invertible to yield the full outflow density distribution p(v,x), they can be used to help 
discriminate among different proposals for the nature of the primary wind and outflow 
acceleration mechanism (e.g. [Masson fc Ghernin (1992"}] , |Ghernin fc Masson (1995}| , agar 
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et al (1997)). Models in which the outflow consists of a shell of ambient matter swept up 
by a wide-angle radial primary wind in a momentum-conserving fashion appear able to 
account for many observed outflow properties, provided appropriate stratification exists in 
both wind and ambient media ( [Shu et al (1991)| , [Li fc Shu (1996)| ). In alternative models, 
the primary wind is assumed to be jetlike (axial flow velocities) and the molecular outflow 
is the manifestation of a bow shock in the ambient medium (e.g. [Blondin, Fryxell, fc Konigl 



T990]i [Masson fc Chernin (1993]| , [Raga fc Cabrit (1993)[ , [Stone fc Norman (1993) |, [Chernin 



et al (1994)1 ). In this case the component of the outflow momentum transverse to the jet is 



driven by pressure gradient forces in the working surfaces at the jet head and within the 
jet (when variability leads to internal shocks), and the agreement with observed outflow 
properties depends very sensitively on the cooling rate of shocked gas. 

The present uncertainty about outflow acceleration mechanisms highlights (and derives 
from) a more basic uncertainty regarding the origin and nature of the primary wind which is 
ultimately responsible for both the observed jets and outflows. Because observed velocities 
of optical jets are comparable to stellar escape speeds, the jet material almost certainly 
originates near the star. Jets appear quite collimated in density down to distances within a 
few tens of AU of the source ([Ray et al (1996)[ ). The observed optical jets are likely only 



the ionized, strongly-emitting inner portions of a broader, neutral wind with lower density 
and lower outflow speed than the jet (for observational evidence, see e. g. [Heathcote et al 



[(1996)1 , [Hartigan (1997)| and references therein) . The presence of separate high and low 
velocity blueshifted components in forbidden line profiles suggests that these winds may 



originate over a range of radii in the circumstellar disks (partigan, Edwards, fc Ghandour 



1995)1) . An open question, however, is whether the apparent collimation of optical jets is 



due to cylindrical density stratification of a primarily radial MHD wind which originates 
near the star ( ^hu et al (1995)] , [Li (1996'F)[ ), or whether there is an MHD wind originating 



in an extended region of the disk (out to ~ 100 AU) which is itself well-coUimated in both 
velocity and density, and which helps to coUimate the observed jet and drive the larger-scale 
outflow (e.g. [Pudritz fc Norman (1986)] , [Konigl fc Ruden (1993)[ ). If observations end up 



demanding the latter, many uncertainties remain in models for MHD disk winds, especially 
in relation to the requirements for producing a collimated (but not recollimated) flow 
without singularities or other unphysical behavior. In the remainder of this Introduction, 
we describe current ideas about MHD winds driven from rapidly-rotating accreting systems 
( "magnetocentrifugal winds"), and discuss some of the difficulties in previous models of 
disk winds. This motivates the present work, which develops models for magnetocentrifugal 
winds that are constrained to become fully collimated at large distances from the source, 
and describes the general properties of such flows. 
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1.2. Magnetocentrifugal Disk Winds 



The basic physics of steady magnetocentrifugal winds has been outlined in numerous 
sources; for a recent pedagogical review, see [Spruit (1996)| . [Heyvaerts fc INorman (1989) 



and Heyvaerts fc Norman (1996) have used analytic arguments to predict how asymptotic 
streamline collimation should develop in winds with varying properties (an analogous 
treatment for the relativistic case is given by [Chiueh, Li, fc Begelman (199r)|) . The mixed 
hyperbolic and elliptic nature of the general governing equations ( [Heinemann fc Olbert 



|( 1978)1) , however, leads to technical difficulties in formulating and finding exact solutions for 
the steady-state problem in terms of appropriate boundary conditions (see e.g. [Bogovalov 
(1994)1) , due to the presence of critical surfaces - projected to curves in the poloidal 
plane - within the flow. Nevertheless, steady-state, fully two-dimensional solutions have 
been obtained for the specific cases of a "split-monopole" poloidal magnetic field interior 
boundary condition ( |!Sakurai (1985)| , [Sakurai (1987)1 ), and for the case of an "X-point" 
poloidal magnetic field interior boundary condition - strongly pinched magnetic fields 
fanning out of the disk near the star (|Shu et al ( 1994b )| , INajita fc Shu (1994)| , |Shu et al| 



1995]]) 



To circumvent the difficulties associated with free critical surfaces, an approach 
that was earlier adopted by [Blandford fc Payne (198T)] (hereafter BP) is to look for 
families of solutions with certain pre-determined symmetries. In BP, the fundamental 
assumption is that all velocities obey the same v oc r~^/^ scaling with spherical radius as 
the Kepler velocity, and that the density obeys p oc r~^/^ (the scaling associated with a 
mass-conservative spherical wind radially flowing in a Kepler potential). The magnetic field 
components then must behave as -B oc r~^/'^ in order to have the corresponding Alfven 
speed v\ = B'^/{4np) oc r~^. By assuming these scalings for the radial dependence of all 
quantities, the governing PDEs of the MHD wind problem are converted to ODEs in angular 
coordinate 6 (or, as many workers have framed the problem, in the scaled height above the 
midplane). The assumed similarity scaling results in any critical surface coinciding with a 
radial line {6 = const.). [Blandford fc Payne (198^ obtained a family of solutions which 
successfully navigated the Alfven transition, directly demonstrating that cold MHD winds 
from rotating disks can reach large velocities and collimate their streamlines toward the 
poles. 

The BP r- self-similar solutions have subsequently been extended and generalized 
by other workers, including relativistic solutions by [Li, Chiueh, fc Begelman (1992)[ 
and [Contopoulos (1994)[ , and nonrelativistic solutions with different basic scalings by 
[(Jontopoulos fc Lovelace (1994)[ (hereafter CL) and [b'erriera (1996)[ In addition, [Sauty[ 
[fc Tsinganos (1994)[ and [Trussoni, Sauty, fc Tsinganos (1996)[ have obtained MHD wind 
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solutions with different basic symmetry (latitudinal, rather than radial, self-similarity) via 
alternative choices for separation of variables, for application to flows from rapidly-rotating 
stars rather than accretion disks. 

Eventually, a theory of magnetocentrifugally driven winds should be able to treat 
time-dependent, non-axisymmetric configurations with a fully self-consistent connection to 
the disk or star where the flow originates, and to the surrounding medium. Time-dependent, 
axisymmetric simulations of outflows from accretion disks have begun to be pursued by, 
e.g. lUchida fc Shibata (1985]| , [Shibata fc Uchida (1986)1 , and [Stone fc Norman (1994)1 



(who include the internal disk dynamics), and by lUstyugova et al (1995)1 , IKoldoba et al 



1(1995)1 , and [Uuyed, Pudritz, fc Stone (1997)1 (who treat the disk as a boundary condition 
for the wind). The simulations produced so far verify that well-collimated MHD outflows 
(sometimes with intriguing time variability) can be generated by differentially-rotating 
disks threaded by a mean magnetic field, and that rapid disk accretion can occur as a 
result of angular momentum removal by the wind. However, the major technical effort 
required to produce these simulations (as well as the steady fully-2D solutions cited above) 
makes it difficult to explore parameter space extensively. Thus, it remains quite useful to 
study the properties of magnetocentrifugally driven outflows with models that impose a 
symmetry in advance to reduce the problem to coupled ODEs. The r-self-similar ansatz 
may be particularly appropriate for modeling outflows from accretion disks with a large 
dynamic range of radii, in which the flow may approach a scale-free solution well away from 
the boundaries. Indeed, the simulations of lOuyed, Pudritz, fc Stone (1997)1 demonstrate 



that steady, r-self-similar winds may naturally develop when a mass flux is driven from the 
surface of a magnetized, rotating disk. 



1.3. Wind-Driven Accretion 

In general, accretion in a disk may be driven in part by local stresses (acting on radial 
scales comparable to the disk thickness, and often parameterized by an "a" viscosity, e.g. 
Pringle (1981)1) , and in part by larger-scale forces such as those associated with global spiral 
density waves (e.g. ISpruit (1987)1 ) an MHD disk wind (e.g. BP, Pudritz fc Norm"anl 



1(1986)1 ). An interesting special case of an accretion disk/wind system is the extreme one 
where the disk wind from a given annulus removes the angular momentum and energy 
needed for what is left of the disk to accrete to the next annulus closer to the star. If we 
assume a steady state for both wind and disk, the conservation equations yield a differential 
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mass transfer rate \dM\Y/dR\ = \dMD/dR\ at radius R, with 

dlnMw f J .\ I 



1 =-r, 1 



d\nR \nR^ J 2' 
in terms of the specific angular momentum J in the wind and the local Kepler rotation 
rate of the disk Q (this expression assumes a thin disk in Kepler rotation). For self-similar 
solutions where density p oc r~'^ in both the wind and disk, this would imply that the wind 
requires a ratio of the cylindrical Alfven radius Ra to the streamline footpoint radius Rq of 

1/2 

I (2) 



Ra /4 - 2q" 



Re \3-2qJ 

because J = flR\ (see §2.1). Solutions which lose very little mass in the outer disk (g nearly 
3/2) must carry large specific angular momentum, and conversely, solutions with massive 
winds from their outer disks (small q) must have relatively low specific angular momenta in 
their winds if they are to be consistent with a steady state inflow/outflow. 

Even if the Alfven radius does not satisfy equation (Q), a self-similar MHD disk wind 
will drive inflow in the disk with the disk accretion rate Mij{R) a power-law in R. For these 
more general cases, the surface density at any point in the disk would either increase or 
decrease in time depending on whether Ra/Ro is greater or less than the value in equation 
(@); the local mass deposition/removal timescale would be greater than the local accretion 
timescale by a factor 2 [3 - 2g - l/((i?A/i?o)^ - 1)]~^ 

The connection between mass and angular momentum loss in the wind and accretion 
in the disk was explored by [Konigl (1989)| and |Wardle fc Konigl (1993)| , using local models 
for the disk and connecting to the BP wind solutions; they showed that the field geometry 
required for a disk wind can be self-consistently provided by a diffusive (ambipolar or 
Ohmic) disk. The work of P!ji (1995)| and [Ferreira fc Pelletier ( 1995)| incorporating additional 
dynamics supported these conclusions. Most recently, [Ferriera (1996)| has obtained steady, 
r-self-similar combined inflow/outflow solutions for the case of ordinary resistivity in 
the disk region, while [Li (1996£i]| has obtained steady inflow/outflow solutions assuming 
ambipolar diffusion in the disk (although neither set of solutions treats the wind far from 
the disk in a complete fashion; see below). In both of the last two cases, the solutions have 
q > 1.4, corresponding to a relatively small fraction of the disk mass being lost to a wind. 



1.4. Critical Points and Wind Asymptotics 



An outflow that accelerates from low velocity near the surface of an accretion disk 
to high velocities at a distance (permitting escape from the gravitational potential) may 
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pass through several points where the flow changes physical character. Such transitions 
occur when the flow speed surpasses the speed of an allowed wave mode (sound waves 
for unmagnetized or MHD waves for magnetized winds), and are manifest by apparent 
singularities in the equations governing the flow. Accounts of the nature of critical points in 
general and restricted MHD winds are given, e. g., by [Heinemann fc Olbert (1978)| , Spruit 



1996)1 , P?singanos et al (1996)| ; we summarize some of the main points here. 



An important issue is whether, and where, the equations pass from being elliptic to 
hyperbolic at large distances from the source. Fully two-dimensional (axisymmetric) MHD 
flows become hyperbolic when the poloidal speed Vp exceeds the fast-mode MHD wave 

speed in the the poloidal direction Vf = (1/2) |cg + vX+ (c^ + v\y — 4Cgf^ p ^ | (where 
Cs is the sound speed, va is the Alfven speed, and f A,p is the Alfven speed associated with 
the poloidal magnetic field component Bp). For a cold flow where Vf — > va, the fast-mode 
Mach number for the poloidal flow is Mp = Vp/vA', the 2D equations become hyperbolic 
when Mp > 1. The full 2D equations also have an apparent singularity at the Alfven 
critical point where fp/fA,p = Ma = 1; here the fluid speed in the poloidal direction equals 
the wave speed of the Alfven mode propagating in the poloidal direction. 



As first pointed out by BP and recently reexamined in detail by p?singanos et al (1996' 



the restriction to an r-self-similar model alters the nature of the critical points. With 
the assumption of r-self-similarity, the poloidal-plane PDE for cross-field force balance is 
converted to a second-order ODE in 6. Critical points in the equations now occur when 
the coefficient of the highest-order ^-derivative in the cross-field equation passes through 
zero. The modified fast-MHD and Alfven points occur in the r-similar MHD equations, 
respectively, if \ve\/vf^ e = 1 (where Vf^ g is obtained by replacing VA,e for va,p in the above 
expression for Vf), and if \vg\/vA,e = 1- For a cold flow, = Vf^ g = va', the fast-MHD 
wave speed is the same for all propagation directions. For an r-similar flow the modified 
fast-MHD critical point would occur at MF(|t>6i|/fp) = \vg\/vA = 1, beyond the point 
Mp = Vp/vA = 1 where the fast-mode transition occurs (there are now no singularities in 
the equations at the location where Mp = 1). Because \vg\/vA,e = \vp\/vA,p, the Alfven 
point in the r-similar flow still occurs at Ma = 1- 

Physically, the modification of critical points happens because r-self-similar model flows 
remain a priori in effective causal contact in the f direction (as well as in the direction 
for assumed axisymmetric flows) for arbitrary Mach number. Thus critical transitions 
respect just the projection of the velocity in the direction of the spatial degree of freedom, 
6, relative to the propagation speeds of the various wave families. In particular, effective 
"ellipticity" in the r-self-similar reduced equations is maintained as long as the 9 component 
of the velocity is smaller than the fast-MHD wave propagation speed in the 9 direction. 
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Vf,e (= "^A for a cold flow). Geometrically, the condition is that the r-similar flow becomes 
effectively "hyperbolic" when the projection of the "minus characteristic" in the 9 direction 
becomes negative (cf. [I'singanos et al (1996)| ; [Contopoulos (19lJF)D . The corresponding 2D 



equations remain formally hyperbolic (in the sense that characteristic curves are deflned) in 
the whole region of the poloidal plane where Mp > 1. Since, however, only the 9 direction 
is a spatial degree of freedom, only the projections of the characteristics in the 9 direction 
are relevant for the propagation of information from the boundaries into the body of the 
solution. Thus only if the flow were to pass through and remain at \ve\lvp^ > 1 could we 
regard the information about the nature of boundary conditions at the pole as irrelevant 
to the solution. Even if such a solution were found, its relevance to real winds would be 
questionable because the propagation of information from the boundaries in the radial 
direction is ignored a priori. 

More generally, full speciflcation of the steady magnetocentrifugal wind problem 
requires both the flow equations and a statement of the desired boundary conditions. 
Whether or not the flow becomes hyperbolic, the interface with the ambient medium inside 
and outside the flrst and last streamlines of the outflow must, at least in principle, govern 
the latitudinal extent of a steady state wind (e. g. |Shu et al (1995)| ). Since, however, the 



r-self-similar model is inflnite radially - with inflnite flux - there is no proper interior or 
exterior to the wind. The only solution boundaries for a cold, r-similar wind lie along ^ = 
and 9 = 7r/2. The boundary at = 7r/2 corresponds to the disk surface, where the solution 
should be able to match to a subsonic flow; several authors have discussed the additional 
constraints this imposes (see §1.3). 

The boundary at = corresponds to the asymptotic limit of the flow at a large 
distance from the source. Even if an r-self-similar solution covers all angles, hence 
formally fllling space and leaving no room for an ambient medium to match, a solution 
should be "physically reasonable" in the sense that if spatially truncated, the solution 
could be embedded within a larger (non-self-similar) solution without any very particular 
requirements. Thus, we would like to obtain self-similar solutions where the ^ — > 
asymptotics could sensibly be matched to an ambient medium with generic properties. Such 
a matching would select among possible input parameters those which yield acceptable 
asymptotic solutions. 

Previous r-similar MHD wind solutions containing both toroidal and poloidal flelds - 
for a variety of scaling parameters - have not explicitly implemented matching to boundary 
conditions at the pole. In some work, numerically-integrated solutions are halted when 
they become singular near the modifled fast-MHD point at \vg\/vA = 1, since no regularity 
condition is applied (e.g. [Ferriera (1996)| ). Other solutions are numerically integrated 
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away from the equator and halted at an arbitrary point near the pole, with jfel/fA < 1 
throughout the computed region (e. g. BP and CL). A difficulty with these solutions is 
that many implicitly require special boundary conditions. Namely, a large class of solutions 
(with q > 1) generically recollimates (i.e. the cylindrical radius i? of a streamline reaches a 
maximum and turns around), and another large class of solutions (with q < 1) generically 
shows radial oscillations of the streamlines (CL; see also |Sauty fc Tsinganos (1994]|) . 



Contopoulos fc Lovelace (1994)| found a third class of solutions (with g > 1) which does 
not show recollimation within the numerically- integrated regime, but they do not explicitly 
implement a ^ asymptotic boundary condition. 



On the other hand, [Contopoulos (199F)| has found a related solution with purely 



toroidal fields (where the flow does not accelerate from zero speed at the disk but instead 
requires a large initial poloidal injection speed, and since f A,p = there is never an Alfven 
transition) with the input parameters tuned such that the solution makes a transition 
through the modified fast-MHD point \ve\/vA = 1- While solutions containing both toroidal 
and poloidal fields can in principle make a transition through \vg\/vA = 1, examples of such 
solutions have not yet been obtained due to the additional numerical effort demanded. The 
added requirement of making the \vg\/vA = 1 transition would lower the number of free 
parameters in the solution by one, so that, for example, the two-parameter family of cold 
wind solutions of BP would become a one-parameter family. 

In the present work, we set out find MHD wind solutions which do match to a specific 
asymptotic prescription. Inspired by observations which show strong coUimation in both 
p and V, we seek self-similar solutions which have cylindrical asymptotics. Thus, while 
our basic set of dynamical equations is equivalent to the cold wind equations of CL or 
Kerriera (1996)| (although mathematically represented quite differently), we supplement 



them with an additional boundary condition at the pole (^ = 0). The solutions we find are 
not global, in the sense that they do not directly connect onto the interstellar medium, 
or to the medium interior to the outflow (e.g. hot plasma, axial fields, or a mass-carrying 
MHD wind originating from the central object/inner disk). Nevertheless, the uniform 
cylindrical asymptotic behavior of these solutions allows a conceptual matching onto simple 
boundary conditions of constant (high) pressure interior, and (low) pressure exterior, to 
the region of the outflow. The levels of the ambient pressure would then select among 
the possible solutions by matching pressures at the boundaries. Although the condition 
of exactly cylindrical asymptotic collimation is probably more extreme than most real 
jet/outflows, the properties of these solutions makes clear the stark contrast in behavior 
between coUimated and uncoUimated magnetized flows. 
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1.5. Plan of Paper 

In §2, we first (§2.1) state the governing equations for a steady, axisymmetric, MHD 
flow, then (§2.2) present an exact analytic solution for self-similar, differentially- rotating 
cylindrical flows which represent the chosen asymptotic boundary conditions far from the 
disk, and finally (§2.3) detail the r-similar wind ansatz and the resulting reduced forms 
of the wind equations. §3 presents solutions for disk winds that satisfy the self-similar 
equations with cylindrical asymptotics, and §4 discusses the properties of these solutions 
and compares to previous models and observations. Various mathematical details are 
described in appendices §§A.1-A.4 



2. Magnetocentrifugal wind equations 

2.1. Steady, axisymmetric ideal MHD equations 

Consider a steady, axisymmetric magnetized flow with negligible resistivity. From 
axisymmetry, the poloidal magnetic field can be expressed in terms of a flux function $ as 

B,^^^. (3) 

where R = r sin 6 is the cylindrical radius. From the z component of the induction equation 
V X (v X B) = 0, it can be shown that $ is conserved on streamlines, Vp • V$ = 0, and that 

Bp = (3pvp (4) 

for some f3 - i.e. the poloidal streamlines and fieldlines are parallel. Combining the 
continuity equation V ■ (pv) = and V ■ B = shows that f5, which represents the ratio of 
magnetic field to mass flux, is also conserved on streamlines: (3 = /?($) and Vp ■ V/? = 0. 
Defining 

n^[v^-Bj{(3p)]/R (5) 

and 

u = V — nR(p, (6) 

we obtain B = jSpu, i.e. the flow in a frame rotating at Q is frozen to the local field line. 
Using V X B = f2V$ and the induction equation, one can show that Vfi x V$ = 0, so 
that Q is constant along field lines. Thus f2($) is the rotation rate of a given field line, and 
on the field line labeled by $ in the frame rotating at Q, the flow and field are parallel, 
u II B. Each field line, of course, can have a different rotation rate fi($). From equation 
(|), Up = Vp; we use the symbols inter changably. 
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The toroidal component of the momentum equation yields an additional conserved 
quantity along field lines, the specific angular momentum J associated with the combined 
matter fiow and Poynting fiux 

J=m = R{v,-^). (7) 
Using the definition we write 

j = viR^ + Ru^ii - Mi^) = nnl, (8) 

where we use the definition for the Alfven Mach number Ma of the flow 

Note that fp/fA,p = u^/vA,ip = |u|/fA = Ma- The Alfven radius R = Ra is the point where 
a flow makes a transition from Ma < 1 near the disk to Ma > 1 in the wind. 

Taking the component of the momentum equation in the direction of B yields the 
Bernoulli equation, which states that 

£ = i|u|2 + + K - -n'^R^ (10) 

is constant on fleld lines, S = £"($). Here, we have assumed an ideal gas with ratio of 
speciflc heats 7, and throughout this work we shall assume the gravitational potential 
Vg = GM/r of a central point mass. In the absence of heating and cooling, the flnal 
conserved quantity on fleld lines is the speciflc entropy, such that K = Pp^'^ = K{^). 

The flnal dynamical equation is obtained by taking the component of the momentum 
equation along V$. The resulting Grad-Shafranov equation describes force balance in the 
direction perpendicular to the poloidal fleld Bp, and can be written 



fl-M^)^ 



B^dJ ^ d£ |BPc//5 dVL dK , , 



In this work, we will flnd it convenient to express the various flow and fleld variables 
in terms of $ and the conserved quantities /?, J, f2, and K. The toroidal speed in the 
rotating frame is obtained from equation (S) as 



JjR -nR 
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the inertial-frame toroidal speed is 

j/R - Mi^nn 

The magnitude of the poloidal speed can be obtained from equation (^) as 

An alternative expression for the poloidal speed comes from the definition (|^) together with 
Up = 'Bp/ (3 p. Equating the two expressions yields an equation for p. For the case of a cold 
flow with P = 0, p is the solution of a quartic equation in terms of $, V$, /3, J, fi, £ and 

r, e. 



2.2. Self-similar magnetized cylindrical flows 

In this section, we consider the steady, ideal, axisymmetric MHD equations without 
gravity, and find a class of solutions representing rotating, cylindrically self-similar flows. 
In these solutions, the poloidal velocity is everywhere axial (vp || 5), all speeds scale 
with cylindrical radius = rsin6' as f oc the density scales as p oc R~'^, and 

all components of the magnetic field scale as -B oc R~^^^'''>/'^ . All flow variables are 
independent of z. The cylindrical solutions have identical similarity scaling in sin 6' to 
their similarity scaling in r. These cylindrical solutions define the asymptotic angular 
behavior for the radially self-similar disk outflows considered in this work, which satisfy the 
general r-self-similar equations to be presented in § p.3| (i. e. the same radial scaling as the 
asymptotic state but arbitrary scaling in 9). For simplicity, we specialize here to cold flows 
(P = 0), but the more general expressions are easily derived. 

Starting from the above ansatz for scalings, any solution must have $ oc R'^ , 
n oc p oc P«/^ J oc Pi/^ and £: oc P-^ The Alfven Mach number Ma is uniform 

throughout the flow. In order to normalize these power-law solutions, we choose a fiducial 
field line $ = $i lying along P = Pi, and rotating at a rate Qi = f2($i). We define 
j = J/VlR^, e = —£/ {VIRY , and m = M^^ = all of which quantities are independent of 
P. Then, normalizing all speeds by VlR = fiiPi(P/Pi)"^/^, we have inertial-frame toroidal 
velocity 



(15) 



rotating-frame toroidal velocity 



_ J - m 
OR ~ 1 -m' 

(16) 



nR~ 1-m'' 
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poloidal speed Vp = Vz 



QR 



1 -2e 



(j - 1)^ 



(1 — m)2 



1/2 



(17) 



and toroidal and poloidal Alfven speeds 



and 



(19) 



For the cold, zero-gravity solutions, the Bernoulli equation (|T0|) with the cylindrical 
ansatz reduces to 

|U|2 



l-2e, 



(20) 



so that v\/{QRy = m{l — 2e). The only remaining constraint on the flow is the 
Grad-Shafranov equation (pUj), which for our present assumption of a cylindrically-stratified 
flow is just the R component of the momentum equation. 



0-^ 



1 djR^Bl) 1 d{Bl) 



R SirpR^ dR Snp dR 



(21) 



where we have dropped the thermal pressure term assuming a cold flow (i.e. Cg small 
compared to the flow and Alfven speeds). Now employing the self-similar scaling for the 
magnetic field |B| oc this reduces to 



n 2 , 9-1 2 1+^2 



(22) 



Substituting in equation (|2^ ) for u^, va,^, and va in terms of e, j, and m, we arrive at the 
quadratic equation 



= m2(l-2e) 



1 + q 



m 



1-q 



;i + g)e 



J 



(23) 



The physical solution has 



.-1 1-9 



+ (1 + g)e + 



1-q 



+ (l + g)e +2f (l + g)(l-2e) 



1/2' 



(24) 

Thus, given q and the values of the angular momentum and Bernoulli parameters j and e, 
we obtain cylindrical solutions by substituting equation ( plj) in equations (|15l)-(|19]). The 
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fast-mode Mach number Mp = |vp|/fA is also constant throughout the cyhndrical flow, and 
is written 



When I — q >> j, e, an approximate solution to equations ( p!5| - |25D is v^/{QR) ^ j, 



vJinR) ^ (2j - 2e)V2, VA,J{nR) ^ j[2/(l - q)]'/^ v^jm) ^ 2j(j - ef'\l - q)-^/\ 
Ma ~ J"M(1 - ?)/2]'/', and Mp ^ j'^l - qf'^j - efl\ 



These self-similar cylindrical solutions only exist for a limited range of parameters. In 
particular, the original Blandford-Payne scaling p oc r~^/^, B oc r~^/^ is not compatible 
with the self-similar cylindrical asymptotic solutions described here. In fact, equation (^) 
can only be satisfied only for g < 1, i.e. when the density and magnetic field profiles are 
less steep than (the same criterion also holds when thermal pressure is included). 
Physically, this is true because the tension associated with the toroidal field is the only 
inward force that can oppose the outward centrifugal force and outward (for g > — 1) force 
associated with the magnetic (and thermal) pressure gradient. Only for g < 1 is the hoop 
stress large enough to enforce cylindrical collimation, and the closer g is to 1, the larger 
the collimation radius relative to the launch point. Notice that it is only the inclusion of a 
toroidal field f a,(^ 7^ that permits cylindrical solutions for fiows in which the magnetic 
pressure decreases outwards. In cases where the magnetic pressure (i.e. the fiuid energy 
density) increases outward (g < —1), both the poloidal and toroidal components of the 
field apply stresses that oppose the centrifugal force in equation (p2D, and collimation 
occurs at relatively small radii. However, we do not consider cases with g < — 1 likely to be 
appropriate models for winds from extended accretion disks, so we do not consider them 
further herein. 

The axial current carried by the self-similar cylindrical fiow scales as R'^ . Therefore, 
for g < 1 the current increases with the radial scale of the fiow, and there is no singularity 
at the axis. The total mass, momentum, and energy per unit time carried by the fiow 
within R scale as M oc i?^/^""^, P oc R^~'^, and E oc i?^/^"''. Thus the cylindrical solutions 
with g < 1 have mass and momentum fiows dominated by the outer regions, and energy 
fiow dominated by the interior (exterior) for g > 1/2 (g < 1/2). 

A further constraint on the parameter space for which self-similar cylindrical solutions 
exist is that the solution of equation ([T7|) be real. The boundary of parameter space is 
found by equating |vp| = and using equation (0); solutions exist for 
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(25) 



2e 



1 + 



4e 1 



(26) 



- 1 + (1 - 2e)i/2 
All solutions must have e < 1/2. 



(1-g) (l + (l-2e)V2)(i_2e)i/2 
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Fig. 1. — Parametric dependence on Rq/Ri and Ra/Rq of the asymptotic self-similar 
cylindrical solution family with q = 0.75. For a given cylindrical solution with field line 
rotation rate fl, specific angular momentum J and Bernoulli constant £ on the streamline 
through Ri, Rq and Ra represent the cylindrical radii of the footpoint and Alfven point a 
disk outflow would need in order to have the same J — QR\, and £ — 

-(3/2)(ni?o)^. (a) Alfven Mach number Ma = Vp/vA,p- Contours show Ma = 1,2,4,8, 
from right to left, (b) Fast MHD Mach number My = Vp/v^- Heavy contours show 
Mf = 0,1,2,4,8,... from right to left; light contours Mp = 0.1,0.2,0.3,.... (c) Poloidal 
speed relative to Kepler speed on the disk. Heavy contours show v^/VLRq = 0, 1,2,3, from 
right to left; light contours show Vp/flRo = 0.1,0.2,0.3,.... (d) Toroidal speed relative to 
Kepler speed on the disk. Contours show v^/VLRq = 0.1, 0.2, 0.3, .., 0.6 from left to right. 
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If we think of a given cylindrical solution as the asymptotic limit of an outflow from 
an accretion disk, we can express the parameters j and e in terms of important physical 
scales in the flow. From equation j = {Ra/Ri)^, where i?A is the Alfven radius for the 
streamline that asymptotes to Ri. From equation (|TU[), the Bernoulli parameter is equal 
to its value where the flow leaves the disk. If the same cold flow streamline originates in 
a Keplerian disk at launch radius Rq, u = on the disk surface so £^ = —{3/2){QRqY . 
Using fl'^ = GM/Rq for the fiducial streamline, e = (3/2)(i?o/-Ri)^. The flow speeds v^p 
and Vz relative to the Kepler speed QRq at the footpoint of any streamline are found by 
multiplying the right-hand-sides of equations (^) and ([TTD by Ri/Rq. From equation (p!^), 
the ratio of the final speed to the Kepler speed at the launch point is bounded by 



< 



2 



(Ra 
\Ro 



1/2 



(27) 



Figure 1 shows the dependence of Ma, Mp, v^p, and Vp on the parameters 
Rq/Ri = (2e/3)^/^ and Ra/Rq = [3j/(2e)]^/^ for an example of the cylindrical solution 
family with q = 0.75. At other values of q, the behavior is similar, with the right-hand 
"corner" of the solution moving toward smaller Rq/Ri and larger Ra/Rq as q increases. 

From equation (p6|), j > e so a cold flow originating in a Kepler-rotating disk must 
have {Ra/RqY > 3/2. Recalling from equation (0) that a fully self-similar inflow/outflow 
solution would have {Ra/ Rq)^ = (4 — 2g)/(3 — 2q), we note that this implies such solutions 
are only possible when q > 1/2. From equation (^), steady inflow/outflow solutions would 

1 /2 

need Vp/{QRo) < [{2q — 1) /(3 — 2g)] ' ; since the maximum value of the right-hand side is 
one (for q = I) the asymptotic axial wind speed for such solutions could not exceed the 
Kepler speed at the launch point QRq. 



2.3. Self-similar steady wind equations and nondimensionalization 



Section ^]T] presents the general ideal MHD equations. Here, we specialize to the case 
of self-similar flows in which all the flow variables are power-laws in the spherical radius 
r, and in particular all velocities are required to follow the Kepler-law behavior v oc r~^/^. 
With this ansatz, we write the functional forms of the density, magnetic flux, field-line 
rotation rate, magnetic-to-mass flux ratio, total specific angular momentum, and Bernoulli 
parameter as 



(28) 
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$ = $1 



' r sin 6 ' 



(3-<7)/2 



(29) 



-3/2 



^(^), 



(30) 



(31) 



J= Ji 



r sin 



1/2 



(32) 



and 



r sm ( 



(33) 



In the above expressions, we have anticipated the asymptotic [6 0) functional form 
of the flow variables by explicitly including the appropriate power-law dependence in sin^ 
in each definition. So that each quantity is normalized by its respective asymptotic value 
on = rsin^ = i?i (for 6 ^ 0, r ^ oo, ), we set n(0) = 1, 0(0) = 1, cu(0) = 1, 6(0) = 1, 
£(0) = 1, and e(0) = 1. Q As introduced in § p^ , we define the constants 



and 



J 


1 PlPi 






-8 


-Si 







(34) 
(35) 



(36) 



for cold flow originating at _Ro in a Kepler-rotating disk, j = {R^/Ri)^, e = (3/2)(i?o/-Ri) 
and m is given by equation (E^. It is also convenient to define the constant 



h = 



{QiRi 



2e 



U - 1)' 



m] 



(37) 



^For a more general radially self-similar flow - not requiring cylindrical asymptotics - we could use the 
same form of the equations and instead normalize by fixing the values of n, 0, uj, etc. on an arbitrary — 9i, 
r — ri. 
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where fp i is the poloidal speed on the streamhne asymptotic to Ri, and the final equahty 
is derived from the asymptotic cyhndrical solution (see eq. |T^). 



From § p.l| , the functions fl, (3, J, and £ are all field line invariants, such that = i^i 
on $ = $1, etc. From expressions (^) - (|33D , this implies that radially self-similar fiows 
have 

UJ = (/)-3/(3-?), b = £ ^ 01/(3-.)^ and e = (/.-2/(3-«) (38) 

for all e. 

To derive an equation for the remaining fiow function n{6), we use the self-similar 
expressions (^) - ( |33D in the equations (0), (^, and ([M]) to obtain 



"Q'W-o/j£)\.,,-2e. + .^-f v. (39) 



(6n)^ \RiJ \1 — mnb'^ ^ 

Here, 

Q2(^) = 0^sin'^+ (^0cos^ + ^^0'sin^ ) , (40) 
where we have set the field line rotation speed QoRo at the footpoint Rq equal to the Kepler 




speed ^GM/Rq, and we drop the thermal energy term. Equation (^) is a quartic equation 
for n{6) in terms of 0(6'), 0'(6') = d(j)/d6, 6, and constant parameters. 

To complete the set of wind equations, we must restate the Grad-Shafranov equation 
(p!T|) in terms of the reduced self-similar functions defined in equations (^) - (|33D. The 
result is a (somewhat complicated) second-order O. D. E. which is linear in (t)"{9)\ the full 
expression is given in § |A.lj of the Appendix. 



As described in § |1.4| , a self-consistent r-similar wind solution must either become 
effectively hyperbolic near the pole (with parameters tuned such that the corresponding 
critical transition through the modified fast MHD point is made smoothly), or else if it 
remains effectively elliptic, the solution must match some physically realistic boundary 
conditions at the pole (again, by appropriately tuning the input parameters). In the present 
work, we satisfy this requirement by matching explicitly to an analytic solution near the 
pole. We choose the i?-similar cylindrical solution introduced in § |2.2| as the = boundary 
condition for our r-self-similar solutions. As described earlier, this limits our choice of the 
scaling parameter to g < 1. Additional requirements on the asymptotic solution arising 



from this choice of polar boundary condition are described in §[A.2 



For the cylindrical asymptotic solutions of this work. Ma and Mp become constant, 
and \v0\/vp = sin 6*, which approaches zero near the pole - hence \vg\/vA — ^ near the pole, 
and solutions remain effectively elliptic. Furthermore, we have found that Mp remains < 1 
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throughout the flow for all the cylindrically-collimated solutions we obtain, so these models 
would be elliptic when considered as solutions of the full two-dimensional equations, as 
well. Thus, the solutions obtained in this work must negotiate only the first critical point, 
requiring a smooth transition through Ma = 1. The conditions that must be met at the 



Alfven point are detailed in § |A.3 



The numerical solutions in this work are obtained by initiating the integrations at 
the pole, and tuning the input parameters until a smooth Alfven transition is obtained. 
For some solutions, we find it convenient to match onto the sub-Alfvenic part of the 
flow with a numerical integration that starts at the equator; the two solutions then 
overlap in the sub-Alfvenic region. Once the input parameters have been found such that 
the flow integration initiated at the pole negotiates the Alfven critical transition, such 
matching solutions require no additional choice of input parameters at the equator to be 
fully-specified. When using the shooting method and starting from the equator, however, 
it is necessary to search for the correct initial value of (f)' at the equator such that a good 
match is obtained. The governing equations become singular at the equator for a cold flow 
(since oo as f — 0); the limiting behavior of the equations and implications for the 



equatorial boundary condition are discussed in §A.4. 



3. Magnetocentrifugal disk wind solutions 

In this section, we present solutions we have obtained for r-similar cold MHD flows that 
meet cylindrical (i?-similar) boundary conditions at the pole. As explained in §§ |2.2| - |0| , 



the solutions we obtain are described by q (the scaling parameter for which p oc r~'^ and 
B (X r~(^+^)/^), and the two parameters j = {Rx/RiY and e = (3/2)(_Ro/-Ri)^, where Rq is 
the radius of a streamline footpoint in the disk, i?A is the cylindrical Alfven radius, and 
Ri is the asymptotic cylindrical radius. For our numerical solutions, we fix q and Rq/Ri 
and then vary R^/Ro until the conditions for a smooth Alfven transition are satisfied (see 



A. 3). We thus obtain a two-parameter family of solutions. 



An example of a solution for q = 0.75 and Ro/Ri = 0.35 is portrayed in Figure 2a-d. 
We find that the solution requires Ra/Rq = 2.412438 in order to make a smooth Alfven 
transition; the Alfven surface lies at 6 = 18° with respect to the pole. From Figures 2a and 
2c, collimation becomes very rapid after the Alfven transition has been made. From Figure 
2c, it is apparent that there is no radial (R) oscillation of the streamlines; instead, in this 
and our other solutions the cylindrical radius R of any streamline secularly increases with 
height above the disk. From Figure 2b, notice that near the disk (at large 6) ua.c/? rises above 
^Ap even before the Alfven transition is made. This effect occurs in all of the solutions 
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Fig. 2. — Flow solution for q=0.75, Rq/Ri = 0.35, Ra/Rq = 2.41. (a) Streamlines (solid 
curves) for equal increments in M. The innermost streamline originates dX R/ Ri = Rq/ Ri = 
0.35 on the abscissa and asymptotes to R/Ri = 1 at infinite Z; the corresponding flow makes 
an Alfven transition dXR/Ri = Ra/ Ri — 2A1. The dashed fine shows the locus of the Alfven 
surface, (b) Fluid speeds (sohd curves) and Alfven speeds (dashed curves) of poloidal and 
toroidal flow/field components, in units of the Kepler speed floRo at the footpoint of a 
streamline, as a function of angle 6 with respect to the pole, (c) Streamline radial distance 
(solid curve) and height above the disk (dashed curve) in units of the footpoint radius Rq 
at the equator, (d) Alfven-mode (solid curve) and fast-mode (dashed curve) Mach numbers 
for the flow. Ma = Vp/vA,p and Mp = Vp/vA- 
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we have obtained. Thus, the acceleration in these flows may be thought of as largely due 
to the gradient in the toroidal field pressure, rather than primarily due to the centrifugal 
force in a rotating, nearly rigid, poloidal magnetosphere (cf. e.g. Spruit 1996). From Figure 
2d, notice that the flow is sub-fast-MHD (Mp < 1) throughout. This is also true for the 
rest of the solutions we have found. Another feature of Figure 2b which holds for the other 
flow solutions as well is that the asymptotic rotation speed exceeds the asymptotic poloidal 
(axial) speed. The slight decrease in poloidal speed Vp near the pole is also a general feature 
of our solutions. 

While the general properties of the flow depicted in Fig. 2 is representative of all the 
solutions we have found, the quantitative characteristics of course differ for each solution. 
Figure 3a shows the relationship among the three parameters q, Rq/Ri, and Ra/Rq that 
characterise flows which become asymptotically cylindrical, for a number of solutions we 
have obtained which successfully make a smooth Alfven transition. Numerically, Ra/Rq is 
obtained as an eigenvalue for any (g, Rq/Ri) pair. From Figure 3a, notice that the solutions 
all hug near the limiting locus {v-p — 0) for cylindrical asymptotic solutions to exist. For 
flxed similarity scaling power g, solutions with relatively more angular momentum (large 
Ra/Rq) collimate relatively close to their footpoints (large Rq/Ri), and also have their 
Alfven radii relatively closer to the asymptotic cylindrical radius (large Ra/Ri)- For flxed 
ratio Ra/Rq of the Alfven cylindrical radius to the streamline footpoint radius, an increase 
in the central concentration of the magnetic flux and density (i.e. an increase in q) implies 
a smaller value of Rq/Ri - i.e. the solution expands to a larger cylindrical radius before 
coUimating. 

In Figure 3b, we show the corresponding value of ^a, the angular distance of the Alfven 
surface from the pole, for the solutions shown in Figure 3a. For flxed q, solutions which 
collimate relatively close to their footpoints {Rq/Ri large) make their Alfven transitions 
relatively near the pole (^a small). For flxed Rq/Ri, the streamlines for different values of 
q follow very nearby paths, and the variation of the fluid speeds with angle 9 is also quite 
similar for the solutions with differing q. On the other hand, the variation of the Alfven 
speeds with 9 differs signiflcantly for solutions with different q and the same Rq/Ri; as q 
increases, the Alfven speed at a given 9 increases, thereby shifting the Alfven point closer 
to the pole. 

The asymptotic flow and Alfven speeds vary with the flow parameters q, Ra/Rq, and 
Rq/Ri- In Figure 4, we plot the values of these asymptotic speeds for the same set of 
solutions as in Figure 3. From Figure 4a, notice that for a given similarity scaling q, the 
asymptotic axial fluid speed increases as the Alfven point moves further from the streamline 
footpoint {Ra/Rq increases). At flxed Ra/Rq, solutions with increasingly concentrated 
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Fig. 3. — (a) Relation among parameters Ra/ Ro, Ro/Ri, and scaling power q, for several 
outflow solutions. Rq, Ra, and Ri are the respective cylindrical radii of the footpoint, Alfven 
transition, and asymptotic lateral expansion of each streamline; q = — 91ogp/91ogr = 
— (1 + 2dlogB /d\ogr). Triangles show solutions with q = 0.5; squares q = 0.75; pentagons 
q = 0.9; hexagons q = 0.99. The dashed curve shows the lower boundary for cylindrical 
asymptotic solutions to exist, i. e. the locus where Vp = for a self-similar cylindrical flow 
of given q (see eq. |2^). (b) Angle of the locus of the Alfven surface with respect to the pole, 
for the solutions in (a) (solid lines connect the solution points of given q). 
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Fig. 4. — Asymptotic fluid and Alfven speeds for full r-similar solutions with various 
parameters, in units of the Kepler speed at the streamline footpoint. Triangles show solutions 
with q — 0.5; squares q — 0.75; pentagons q — 0.9; hexagons q — 0.99. Solid hues connect the 
solution points, (a) Poloidal (axial) speed, (b) Poloidal Alfven speed, (c) Toroidal speed, 
(d) Toroidal Alfven speed. 
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magnetic flux (increasing q) have lower asymptotic f^. From Figures 4b and 4c, the same 
trends hold for the poloidal Alfven speed and the toroidal flow speed. Figure 4d shows that 
the asymptotic toroidal Alfven speed increases with R^/Ro for fixed q. At fixed Ra/Rq, 
the toroidal Alfven speed increases very slightly with the degree of flux concentration 
(increasing q). Note that none of the solutions displayed here have Ra/Rq sufficiently small 
to correspond to the case of a steady, self-similar inflow driven by the outflow at given q 
(cf. eq. H). While we have identified some such solutions, the extremely low asymptotic 
poloidal speeds that they yield makes them unlikely candidates to model real outflows. 

Figure 5a,b shows how the Alfven and fast-MHD Mach numbers Ma and Mp for 
the asymptotic limit of the r-similar flow varies with the solution parameters. For fixed 
similarity scaling q, Ma decreases with increasing Ra/Rq, while Mp slightly increases with 
Ra/Rq- At fixed Ra/Rq, an increase in the concentration of the flow (larger q) implies a 
larger asymptotic Ma, but a smaller asymptotic Mp. On the other hand, for fixed ratio 
of initial to final radius of a streamline Rq/Ri (not shown), an increase in q implies both 
smaller Mp and Ma, asymptotically. 

The flow solutions we have obtained are in general magnetically dominated in the 
asymptotic regime, in the sense that for most of the solutions, the energy, momentum, 
and angular momentum fluxes of carried by Maxwell stresses exceed the respective kinetic 
fluxes. These flux ratios are given, respectively, by 





Iv^Ma 


Em 


2QR\vA,,fi\ 



(41) 



^ = i^M,-^ - Mfr, (42) 

and 

^ v^Ma ^^2^ 
Jm \'^a,<p\ 

(Note that the small electric contribution (1/2) M^'^{^lR/cY to the denominator of equation 
(1^) has been dropped.) Figure 6 shows the asymptotic values of these flux ratios for 
the same solutions as discussed above. A few cases with small q and Ra/Rq have kinetic 
angular momentum flux exceeding the magnetic angular momentum flux, but otherwise the 
solutions are magnetically dominated. 
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Fig. 5. — Asymptotic Alfven and fast MHD Mach numbers Ma = Vp/v\^p (a) and 
Mp = Vp/vA (b) for r-similar solutions with various parameters. Triangles show solutions 
with q = 0.5; squares q = 0.75; pentagons q = 0.9; hexagons q = 0.99. Solid lines connect 
the solution points. 
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Fig. 6. — Asymptotic ratios of kinetic to magnetic energy flux (a), momentum flux (b), and 
angular momentum flux (c) for r-similar solutions with various parameters. Triangles show 
solutions with q = 0.5; squares q = 0.75; pentagons q = 0.9; hexagons q = 0.99. Solid lines 
connect the solution points. See text for deflnitions of fluxes. 
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4. Summciry and Discussion 

The calculations presented in this paper explore the proposal that observed narrow 
optical jets from PMS stars faithfully represent an inherent strong collimation in both 
density and velocity of the winds which are magnetocentrifugally driven from these young 
star/disk systems. To this end, we have developed models of MHD disk winds in which all 
streamlines asymptote to lie on the surfaces of a series of nested cylinders, at large distances 
from the disk. To make these calculations tractable, we work within the framework of 
radially self-similar flows, for which the density p, all magnetic field strengths B, and all 
velocities v have power-law scaling in the spherical radius r: p oc r~'^, B oc r^*^^"'"'^^/^, and 
V oc r~^/^, for arbitrary q. The shape of wind streamlines and the dependence of all flow 
quantities on 6 are then calculated self-consistently from the steady-state MHD equations. 
The r-self-similar ansatz is a convenient mathematical idealization, and also may serve as 
a good characterization of the possible flows well away from the boundaries of a disk with 
large dynamic range. Along a ray of ^ = const., the decrease of the dynamical variables (v, 
B, p) with r occurs because the collimation of the flow brings streamlines from successively 
more distant regions of the disk to cross the ray. Since the farther regions of the disk are 
characterized by lower velocity, and by assumption lower density p and magnetic field 
strength B. the wind will reflect these decreases correspondingly. Note in particular that 
the present mathematical formalism allows for fast flows at large distance r from the source 
even though v oc r~^/^: if Vp increases rapidly with decreasing 9 (toward the i-axis), then a 
"jetlike" flow can exist along the pole. 

We begin our analysis, in §2.2, by presenting an exact, analytic family of solutions for 
rotating, cylindrically-symmetric, axial {vji, B^ = 0) MHD flows which are self-conflned by 
gradients in the toroidal magnetic field. All fiow variables are independent of the distance 
along the flow axis z, and have power-law dependence p oc B^ oc Bz (x 

oc Vz (X R"^^^ on the cylindrical radius R — rsinO. The three-parameter family of 
solutions can be described by q together with j and e, where the latter two respectively are 
the speciflc angular momentum, and the Bernoulli constant (equivalent to the fluid energy 
in a rotating frame), in units of local rotation rate of a magnetic fleld line Q and cylindrical 
radius R. These parameters are natural to use when connecting to a flow from a disk in 
Keplerian rotation. This fully-analytic family of shearing, rotating, axial winds includes 
both sub- and super- fast MHD flows, i.e. f ^ < or > va- The presence of a toroidal fleld 
VA,<f 7^ permits conflnement of rotating flows in which the magnetic pressure decreases to 
larger R - i.e., flows in which the outward centrifugal force and outward pressure gradient 
forces are balanced by an inward magnetic tension ("hoop stress"). In order for this to 
work, however, the gradient in the magnetic fleld cannot be too steep: for self-similar flows, 
cylindrical collimation is only possible if B drops with R no faster than R"^, i.e. ^ < 1. 
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We next, in §3, take the family of rotating, shearing cyhndrical flows as the set 
of desired asymptotic boundary conditions far from the disk, and ask what sort of 
magnetocentrifugally-driven disk winds can achieve this fully collimated state. We present 
a range of solutions for various q, e = (3/2)(i?o/-Ri)^ and j = {Ra/Ri)"^, where Rq is the 
position of a streamline's footpoint in the disk, R\ is the Alfven radius, and Ri is the 
asymptotic radius. Note that the self-similarity of the flow makes e and j the same for 
all streamlines. For each q and e there is a unique j for which the flow passes smoothly 
through the Alfven point; thus the family of r-similar disk winds which become cylindrically 
collimated is labeled by two parameters. The asymptotic characteristics shared by all the 
self-collimating wind solutions presented are (1) v^/^Rq ~ few xO.l, (2) v^/QRq ~ few - 
several xO.l, (3) va^^/va^p » 1, and (4) v^/va < 1- Here, Q is the Kepler speed in the 
disk at the streamline's footpoint Rq. In particular, note that unlike in most previous MHD 
disk wind solutions (which, however, have unconstrained asymptotic states), these flows 
do not accelerate to final axial speeds comparable to or larger than the escape speed from 
the potential well where they originated. Thus, for example, if the innermost streamline 
originates near the surface of a PMS solar-type star, then based on the present solutions 
the maximum outflow speed along the central axis of the "jet" would only be a few to 
several tens of kms~^, much lower than the optical Herbig-Haro jets observed to reach a 
few hundred kms~^ speeds. Intuitively, the result (4) may be understood as follows: a 
cold MHD flow must leave the disk at an angle at least 30° from the vertical. In order to 
be fully redirected upwards, VA,ip must be ^ to accomplish the refocusing. But after the 
flow has been redirected, it is left in sub-fast-MHD state. 

Because the fully-coUimated disk wind solutions we have found have relatively low final 
speeds, we believe that observed optical jets from YSOs are unlikely to be manifestations 
of these or similar flows. Among the many possible alternative models, we list a few: 

(1) The wind responsible for the observed jet originates in the disk, and its final velocity 
collimation is large but less extreme than the present solutions. [Contopoulos fc Lovelace| 



(1994)1 have identified solutions with q slightly greater than one which continue expanding 



laterally (in R) up to the limit of their integration in z, and have v^/VLRq > 1. It would 
be interesting to compare the density and momentum distributions in such solutions with 
observed jets. 

(2) A well-collimated MHD disk wind like those computed in this paper comprises an 
unseen neutral wind which surrounds the optical jet and enforces its observed collimation; 
the wind would also help drive outflows by directly sweeping up ambient material far from 
the jet. The jet itself must originate inside of the coUimating wind, and have a more 
efficient acceleration mechanism (perhaps involving relatively strong poloidal magnetic 
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fields to achieve the magnetic "propeller" effect). Since in this case the jet is not required 
to self-collimate its velocity vectors, the flow which produces the jet could emerge from near 
the star in a relatively isotropic fashion. 



(3) As in the model of Shu et al (1994a,b), [Najita fc Shu (1994)| , and (Shu et al (1995]l the 



whole of the wind originates in a narrow region of the disk near the star and flows outward 
with nearly radial poloidal streamlines and speeds comparable to the stellar escape speed 
for the whole flow. At large distances from the star, the velocity field is still nearly radial 
with only logarithmic streamline collimation due to toroidal field stresses. In this case, 
however, the density appears much more cylindrically coUimated than the velocity field, 
and the optical jet may represent just the densest part of the wind near the axis. 

Since the models presented herein are not "global" in the sense of physically matching 
to interior and especially exterior ambient conditions, it is more difficult to assess the 
possible relationship between such winds and observed molecular outflows. A real disk has 
a real inner and outer edge, and thus there must be a first and last wind streamline. If 
we arbitrarily truncate the source of the self-similar wind at the inside and outside disk 
edges -Ro = R* and Rq = Rb, then the first and last streamlines would exit the disk at a 
finite angle[| and asymptote to helically wrap around the surfaces of cylinders of inner and 
outer radii -Rc = {Ri/Ro)R* and i?w = {Ri/ Ro)Rb, respectively. If we assume a central 
star of mass and total mass- loss rate Mw in the wind between and Rd, then we 
can introduce dimensional scales for the density, magnetic field, and velocity to write the 
asymptotic profile in the wind at cylindrical distance R from the axis as 

P = Pw 7^ , (44) 



R 



w 



o X -(l+g)/2 

B = Bw(^ , (45) 



R 



■w 



and 



Here 



v = vw(^) . (46) 



R 



^ \io-mQ jr-^ J \Mq J viooAu; \rJ WrqJ ' 

(47) 



^Note that to maintain force balance near the outer disk edge there would in fact need to be additional 
wind streamlines filling the region toward the equatorial plane. 
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is the density on the outermost wind streamhne, 



B 



w 



lO-^Moyr-i; \MqJ VlOO AU/ \RJ \nRo 
gives the magnetic field components on the outermost streamhne, and 



gives the components of the wind velocity on the outermost streamline; 

is a normalization constant. The values of Rq/Ri, and asymptotic Vz/^Rq, v^/QRq, 
va,z/^Ro, and va^^/QRq are shown for various models in figures 3d and 4a-d of §3. 

The pressure in the wind is dominated by B'^/Stt. While the very innermost part of the 
parent cloud core may have comparable pressure (at density ~ 10^ cm~'^), the outer cloud 
core ambient medium will generally have pressure far below that of the wind. Thus, either 
the wind as a whole or its surface layers would have to expand laterally until the pressure 
matches the ambient medium (which event occurs depends in part on the stability of the 
wind, a question under current investigation). 

If the wind expands as a whole, keeping the mass load and poloidal speed on streamlines 
unchanged, then in the limit of negligible ambient pressure the wind would fill all An 
steradians. The far-asymptotic solution can be computed similarly to the calculation 



outlined in ^hu et al (1995)| , if we assume that the fiow streamlines adjust their latitudes 
6 at each r so that the force associated with the hoop stress (i.e. the gradient in {RB^Y) 
is zero. The resulting streamlines and density contours have qualitatively similar behavior 



to the [Shu et al (1995)| results for a wind originating from a narrow region in the disk 



near the star. In particular, the distribution of wind momentum fiux with angle has 
dP/dcosO ~ C(r)i?Q^~''^^^/(sin6')^ for C{r) a slowly- varying function, and Rq the footpoint 
radius of the streamline that passes through 6 at that r. For q near unity, this is nearly 
the same momentum distribution as that used by |Li fc Shu (1996)| in their calculations of 



the lobe shapes and line profiles for swept-up molecular shells, so their results would carry 



over to the case of "fully expanded" disk winds. As |Li fc Shu (1996)| showed, a momentum 
distribution near dP/d cos 9 oc l/(sin^)^ and toroidal mass distribution Pamb sin^ 6 in the 
ambient medium well reproduces the parabolic outfiow shell shapes and line profile wings 
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dMsh/dVsh oc v~^'^ characteristic of many observed outflows. The coincidence in momentum 
distribution between the present case and the X-wind models occurs because here, the lower 
poloidal speeds on the outer disk streamlines are offset by a relatively greater mass load, 
compared to the X-wind models. Because winds from the outer portions of disks would have 
lower speeds but higher densities than the corresponding X-wind with the same momentum 
flux, the two cases could potentially be distinguished by the emission properties of the 
region where the wind impacts the ambient cloud. For example, finding H2 emission at low 
latitudes in outflow shells would argue against disk-wind models for sweeping up the shell, 
because a shocked low-speed wind would not reach the few-thousand Kelvin temperature 
required. 

If instead of the wind expanding as a whole, we assume that only the surface layers 
expand to match the pressure in the ambient medium, then much of the wind could retain 
its axial poloidal velocity field. The distribution of axial momentum flux in this interior 
portion would obey dP/dR oc R~'^. If this inner wind sweeps up a shell from a surrounding 
cloud core with initial density distribution pamb oc then the distribution of mass with 
velocity Vs in the shell would obey dMsh/dVgh oc 'uj^/^^~^'^^~^"-^^~''^/^^'^''\ For q near unity, or a 
near 2 (as for the |Shu (1977]| singular isothermal sphere), the power in the distribution is 



near —2, in good agreement with the —1.8 dependence cited by [Masson fc Chernin (1992) 



for the high-velocity part of observed outflow line profiles. Without a model for how the 
wind on the outer streamlines expands and interacts with the ambient gas, however, it 
is impossible to tell whether such a "partially expanded" coUimated disk wind would be 
consistent with observed outflow shapes and the distribution of low-velocity gas. 

In summary, we believe that highly-coUimated fast flows as defined by such famous jets 
as IIII30, IIII34, or HHlll are unlikely to represent the whole story of the primary winds 
from PMS star/disk systems. The models of this paper, though by no means an exhaustive 
survey of jet production possibilities, show that it may be difficult to generate strong 
velocity coUimation and fast outflow self-consistently in a magnetocentrifugally-launched 
wind from an extended region of an accretion disk. The narrow appearance of fast optical 
jets from young stars therefore argues in favor of a surrounding wind, either from the whole 
of the disk or from near the star, which sustains the observed jet coUimation and helps to 
drive molecular outflows on large scales. 
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A. Appendix 



A.l. The Grad-Shafranov equation for r-self-similar flows 

The Grad-Shafranov (GS) equation describes force balance in the direction 
perpendicular to the poloidal field lines. With the adoption of a radially self-similar form 
for the fiow (eqs. ^3|), the GS equation becomes a second-order ODE for the reduced 
fiux function 
written 



Denoting derivatives with respect to 9 as primes, this equation can be 

(Al) 



n 
1 



where 



I mnh'^ 2 ^ 



mnb'^ — 1) Imrafo^ 

" / \ 2 
/ ^A,p \ 



Ri 



2ee 



h 



)sin^ 
nb 



GM/R 



(A2) 



and 7^ = 7^l + 7^2 + 7^3. Here, 
3-q V{hnf 



7^l 



4 (j)h sin^ 9 
(q + l)mnb'^ 



X 



2ee 



RiJ 1 — mnb'^ 



,(A3) 



sin 6* 3 — g 



(A4) 



and 



T^q = (1 — mnb 



h 



, ^ 3 — g cos 9 
2 sin^ 



+ 



X 



2{nby 

Rq 

R'l 



/ 2q ' 

■ sin ^ 50 cos H (f)' sin 6' 

V 3-g , 



, , 3 — g cos 9 



2 sin 6* / 3 — g 



3 — g (3 — g)(l — mn62)2 



(A5) 



with 



V = mnb'' 



R 

2 ( ) sin^ + cj^-2ee 



{nby 
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{mnh'^ — 1] 



— I sm e' + a; 



2ee 



+ 



UJ] 



Rj GM/R ■ 



(A6) 



The function Q, which is proportional to the poloidal field strength, is defined in equation 
(^). The constants e, j, m, and h are defined in § |2.3| . The first two are the free parameters 
which must be chosen before initiating an integration, while m and h may be expressed in 
terms of e and j using the asymptotic cylindrical solution (see eqs. El and E7\). 



A. 2. Asymptotic boundary conditions 

The equation ( [AID contains terms that scale as (sin^)~^ and (sin6')~^, and this singular 
behavior near the pole demands special treatment. The method is to approximate 



= 1 + (j)'{o)e + (f)"{o)ey2 (A7) 

near the pole, and then to expand the governing equation (^) in 6 and collect terms at 



each order. The only part of equation (|A1|) that contains terms of order (sin 6') is TZi 



and equation (^) shows that the expression in curly braces in equation (|A3|) contains no 



order-unity terms. This verifies that the cylindrical solution presented in §^]2| is a valid 
limiting solution of the full self-similar wind equations, to lowest asymptotic order. Next, 
we must collect the order (sin^)~^ terms in the expansion of equation ( |A1| ). By requiring 
that the corresponding coefficient be zero, after some algebra we obtain a linear equation 
which we solve for 0'(O) in terms of q, e, and j. Finally, the expansion of ( [AID order 1 is 
linear in 0"(O); rather than performing the expansion analytically, however, we numerically 
solve for the value of 4>"{0), once q, e, j, and 0'(O) are set. After solving for 0"(O), we choose 
a value for 6i (typically 10~^ is adequate) at which to initiate the integration, then evaluate 
(f) at 6i via equation ([A7| ) , set 0'(^i) = 0'(O) + (j)"{0)6i, and begin the numerical integration. 



A. 3. Alfven transition constraints 

As discussed by BP and other authors, equation ([Al|) can become singular if £ = 0, 
which from equation (|A2|) occurs if fA,p = Vp (i. e. Ma = 1) or va = ve- For waves 
propagating along the 6 direction (i.e. along the spatial direction corresponding to the only 
degree of freedom for a radially self-similar fiow), the first case corresponds to an Alfven 
transition in vg, and the second case to a fast MHD transition in vg. To avoid singularities. 
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physical solutions must have TZ = wherever £ = 0; such points become critical points of 
the flow, and the requirement that the flow passes smoothly through these critical points 
constrains the possible values of the parameters entering the solution. It turns out that for 
the solutions presented in this work, \ve\ < va throughout the flow, so there is no fast MHD 
critical point. The Alfven critical point, however, does impose constraints on the possible 
values of the parameters j and e. 



At the Alfven point where Ma = 1, we use the streamline-invariant equations (|3^) to 
evaluate u, b, £, and e in terms of 0(6'a) = 0a- Then, from equation (0), a smooth solution 
requires ji — u!\A = so that 0a = must obtain. Using the definition = mnb"^, 

the reduced density at the Alfven point must be given by ua = j'^^'^/m. The Bernoulli 
equation at the Alfven point becomes 



de 



4u;0' 



A (3 - g)0pV2 



(A8) 



where T) is defined in equation (|A6|) . Finally, the Bernoulli and GS equations at the Alfven 
point can be combined to yield the requirement that 0^ must be a solution of 



(4^) 



2ee 



2i^\^ sine 



RiJ 



+ 



2tp' sin e cos e 
3-9 



L^2-2ee + 2(|^) sin^ 



+ 3a;2 - 2ee 



(A9) 



Thus, the values at the Alfven point of 0, u;, 6, £, e and n are determined from the 
initial conditions at the pole, while the required value of 0a is given implicitly through 
equation ( |A9|) in terms of these known parameters and the unknown position of the Alfven 
point 9a- For arbitrary values of the input parameters e and j, however, 0' at the point 9a 
where = 0a will not satisfy equation (|7^). Thus, one of the input parameters becomes 
an eigenvalue which is determined by the condition that equation ( [A9| ) is indeed satisfied, 
so that the solution is regular at the Alfven point. We evaluate the eigenvalue using the 
shooting method. Specifically, we choose values of q and Rq/Ri (i-e. e = {3/2){Ro/ RiY), 
and guess a value of Ra/Ri (i.e. j = {Ra/Ri)"^)- We then numerically integrate from the 
pole up to the Alfven point, defined by the criterion = 0a. We use the error in equation 
( |A9| ) near the Alfven point to correct the guess for Ra/Ri, and repeat until a satisfactorily 
converged solution is obtained for the super-Alfvenic region (typical error is ~ 10^^ — 10~^ 
at the Alfven point). We can then step over the Alfven point, and continue the integration 
on the sub-Alfvenic side. 
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A. 4. Equatorial boundary conditions 



For integrations which start at the equator 



7r/2), divergence of n{9) (since 



f p — ^ in the cold-flow limit) means we must analytically expand the relevant equations 
in 6 = n/2 — 9. For a given value of Rq/Ri (i.e. e), we have (f){7i/2) = 0e = (Ro/Ri) ^ 
from equation (^9|) , and from equation (|38D we obtain the equatorial values of the other 
streamline invariants. Given an arbitrary value for 0'(7r/2) = 0^, the lowest-order expansion 
of the Bernoulli equation (|39|) yields the approximation for the reduced density function 



n{6) = (mblS)-' {hm^ 



-1 + 3 




+ 



Ra 
Ro 



1/2 



X 




-1/2 



(AlO) 



where here and in the following the "e" subscript denotes evaluation at 6 = 7t/2. From 
this equation, it is clear that solutions require (3 — q)/{2^/3) < |0g|/0e; this condition is 
equivalent to the well-known requirement that streamlines leave the disk at an angle > 30° 
with respect to the vertical (BP). 



Next, given some choice for 0^, we can expand the Grad-Shafranov equation ( [AID 
to lowest order in 6. The result is linear equation which yields 0" in terms of q, e, j, 
and 0g. For fixed q, we can identify a pair {j, e) which satisfies the requirements for a 



smooth transition at the Alfven point when approached from the pole, as described in § |A.3 
We can then iteratively solve for the value of 0g which allows an equatorially-initiated 
sub-Alfvenic solution to match smoothly onto the pole-initiated super- Alfvenic solution at 
the Alfven point. We generally proceed by choosing a value for Rq/Ri and then searching 
for the unique values of Ra/Ri and cf)'^ which yield a smooth match at the Alfven point. 
When sufficiently accurate values are determined, the pole-initiated and equator-initiated 
solutions can generally cross the Alfven point and overlap the solution on the far side, with 
very small errors. 
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